home *** CD-ROM | disk | FTP | other *** search
/ Collection of Tools & Utilities / Collection of Tools and Utilities.iso / edit / jed096_1.zip / SLANG / SRC / SLMATRIX.C < prev    next >
C/C++ Source or Header  |  1994-04-26  |  4KB  |  174 lines
  1. /* Floating point matrix manipulation routines for S-Lang */
  2. /* 
  3.  * Copyright (c) 1992, 1994 John E. Davis 
  4.  * All rights reserved.
  5.  *
  6.  * Permission is hereby granted, without written agreement and without
  7.  * license or royalty fees, to use, copy, and distribute this
  8.  * software and its documentation for any purpose, provided that the
  9.  * above copyright notice and the following two paragraphs appear in
  10.  * all copies of this software.
  11.  *
  12.  * IN NO EVENT SHALL JOHN E. DAVIS BE LIABLE TO ANY PARTY FOR DIRECT,
  13.  * INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES ARISING OUT OF
  14.  * THE USE OF THIS SOFTWARE AND ITS DOCUMENTATION, EVEN IF JOHN E. DAVIS
  15.  * HAS BEEN ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
  16.  *
  17.  * JOHN E. DAVIS SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT
  18.  * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
  19.  * PARTICULAR PURPOSE.  THE SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS"
  20.  * BASIS, AND JOHN E. DAVIS HAS NO OBLIGATION TO PROVIDE MAINTENANCE,
  21.  * SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS.
  22.  */
  23. #include <stdio.h>
  24. #ifndef FLOAT_TYPE
  25. #define FLOAT_TYPE
  26. #endif
  27. #include "slang.h"
  28. #include "_slang.h"
  29. #include "slarray.h"
  30.  
  31. static SLArray_Type *SLpop_float_matrix (void)
  32. {
  33.    SLArray_Type *a;
  34.    int sa = 0;
  35.    
  36.    if (NULL == (a = SLang_pop_array (&sa))) return NULL;
  37.    
  38.    if (a->type == FLOAT_TYPE) return a;
  39.    SLang_Error = TYPE_MISMATCH;
  40.    return NULL;
  41. }
  42.  
  43.    
  44. /* multiply 2 matrices (assumed float) producing a third */
  45. static void SLmatrix_multiply (void)
  46. {
  47.    SLArray_Type *a, *b, *c;
  48.    FLOAT *aa, *bb, *cc, sum;
  49.    int c_handle, dim;
  50.    unsigned int imax, jmax, kmax;
  51.    unsigned int i, j, k, ofs_a, ofs_b;
  52.    
  53.    if ((NULL == (b = SLpop_float_matrix ()))
  54.        || (NULL == (a = SLpop_float_matrix ())))
  55.      return;
  56.      
  57.    /* Now is a is n*m, then b must be m*x.  Result it n*x */
  58.    
  59.    imax = a->x;
  60.    jmax = a->y;
  61.    
  62.    if ((b->x != jmax) || (a->dim > 2) || (b->dim > 2))
  63.      {
  64.     SLang_Error = TYPE_MISMATCH;
  65.     return;
  66.      }
  67.    
  68.    kmax = b->y;
  69.    
  70.    /* Now result will be  imax by kmax 2d array */
  71.    if (kmax == 1) dim = 1; else dim = 2;
  72.    
  73.    if (-1 == (c_handle = SLcreate_array(NULL, dim, imax, kmax, 1, 
  74.                        'f', 0)))
  75.      {
  76.     SLang_doerror("Unable to create array.");
  77.     return;
  78.      }
  79.    c = SLarray_from_handle (c_handle);
  80.    /* Finally!! */
  81.    cc = (FLOAT *) c->ptr;
  82.    bb = (FLOAT *) b->ptr;
  83.    aa = (FLOAT *) a->ptr;
  84.    
  85.    for (i = 0; i < imax; i++)
  86.      {
  87.     for (k = 0; k < kmax; k++)
  88.       {
  89.          sum = 0.0;
  90.          ofs_a = i * jmax;
  91.          ofs_b = k;
  92.          for (j = 0; j < jmax; j++)
  93.            {
  94.           sum += *(aa + ofs_a) * *(bb + ofs_b);
  95.           ofs_a++;
  96.           ofs_b += kmax;
  97.            }
  98.          
  99.          /* cc[i][k] */
  100.          *(cc + (int) (i * kmax + k)) = sum;
  101.       }
  102.      }
  103.    
  104.    SLpush_array (c_handle);
  105. }
  106.  
  107. static void SLmatrix_addition (void)
  108. {
  109.    SLArray_Type *a, *b, *c;
  110.    FLOAT *aa, *bb, *cc;
  111.    int c_handle;
  112.    unsigned int imax, jmax, kmax, jmaxkmax;
  113.    unsigned int i, j, k, ofs;
  114.    
  115.    if ((NULL == (b = SLpop_float_matrix ()))
  116.        || (NULL == (a = SLpop_float_matrix ())))
  117.      return;
  118.      
  119.    /* for the addition to make sence, they must be same type. */
  120.    imax = a->x; jmax = a->y; kmax = a->z;
  121.    
  122.    if ((b->dim != a->dim) || (b->x != imax) 
  123.        || (b->y != jmax) || (b->z != kmax))
  124.      {
  125.     SLang_Error = TYPE_MISMATCH;
  126.     return;
  127.      }
  128.    
  129.    if (-1 == (c_handle = SLcreate_array(NULL, a->dim, imax, jmax,
  130.                        kmax, 'f', 0)))
  131.      {
  132.     SLang_doerror("Unable to create array.");
  133.     return;
  134.      }
  135.    
  136.    c = SLarray_from_handle (c_handle);
  137.    /* Finally!! */
  138.    cc = (FLOAT *) c->ptr;
  139.    bb = (FLOAT *) b->ptr;
  140.    aa = (FLOAT *) a->ptr;
  141.    
  142.    
  143.    /* Probably more efficent if we work in this order */
  144.    jmaxkmax = jmax * kmax;
  145.    for (k = 0; k < kmax; k++)
  146.      {
  147.     for (j = 0; j < jmax; j++)
  148.       {
  149.          ofs = j * kmax + k;
  150.          for (i = 0; i < imax; i++)
  151.            {
  152.           *(cc + ofs) = *(aa + ofs) + *(bb + ofs);
  153.           ofs += jmaxkmax;
  154.            }
  155.       }
  156.      }
  157.    
  158.    SLpush_array (c_handle);
  159. }
  160.  
  161. static SLang_Name_Type slmatrix_table[] =
  162. {
  163.    MAKE_INTRINSIC(".matrix_multiply", SLmatrix_multiply, VOID_TYPE, 0),
  164.    MAKE_INTRINSIC(".matrix_add", SLmatrix_addition, VOID_TYPE, 0),
  165.    SLANG_END_TABLE
  166. };
  167.  
  168. int init_SLmatrix()
  169. {
  170.    return SLang_add_table(slmatrix_table, "_Matrix");
  171. }
  172.  
  173.  
  174.